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ABSTRACT 

The form of the solar meridional circulation is a very important ingredient for mean field flux 
transport dynamo models. Yet a shroud of mystery still surrounds this large-scale flow, given 
that its measurement using current helioseismic techniques is challenging. In this work we 
use results from 3D global simulations of solar convection to infer the dynamical behavior of 
the established meridional circulation. We make a direct comparison between the meridional 
circulation that arises in these simulations and the latest observations. Based on our results 
we argue that there should be an equatorward flow at the base of the convection zone at 
mid latitudes, below the current maximum depth helioseismic measures can probe (0.75 Rq). 

We also provide physical arguments to justify this behaviour. The simulations indicate that 
the meridional circulation undergoes substantial changes in morphology as the magnetic cycle 
unfolds. We close by discussing the importance of these dynamical changes for current methods 
of observation that involve long averaging periods of helioseismic data. Also noteworthy is the 
fact that these topological changes indicate a rich interaction between magnetic fields and 
plasma flows, which challenges the ubiquitous kinematic approach used in the vast majority of 
mean field dynamo simulations. 
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1. INTRODUCTION 

The pressing need for a better understanding of 
space weather phenomenology and its implications 
has provided impetus for ongoing improvements in 
solar dynamo models. Nonetheless a fully detailed 
and consensual picture adequately explaining the ori¬ 
gin and evolution of the large-scale solar magnetic 
field has not yet materialized. Amongst the sev¬ 
eral t ypes o f ext ant mean-field models, (see review 
bv iCharbonneaul [2010 ). a specific class of dynamo 
models known as mean-field Babcock-Leighton Flux 
Transport Dynamos (FTD) has been particularly suc¬ 
cessful in explaining some of the main observa tional 
f eatures of t he s ola r cy cl e, e.g . dWang et all 119911 
DikDaH&_Charbonneau Il999t iNandv fc Choudhuril 
120021 lKa.ra.k~ et a/.l I20T4I1 . These models rely heav¬ 
ily on the advective role of the meridional circulation 
(MC), a relatively weak, large-scale plasma flow that 
its thought to act as a conveyor belt in the convec¬ 
tion zone (CZ), transporting magnetic flux in radius 
and latitude. In these models (running mainly in 
the kinematic regime), the MC is an important de- 
terminant of the pe riod and amplitude of the cycle 
(|Charbonneaull2010l and references therein). 
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From an observational point of view, the merid¬ 
ional circulation is known with a good degree of con¬ 
fidence at low to mid latitudes in the photospheric 
layers of the Sun, through Do ppler measurements and 
track ing of surface features (jlJlrichl 120101 : lHathawavl 
2012). This surface component of the MC appears 
to play an important role in the inversion of the so¬ 
lar dipole. Observations and surface advection mod¬ 
els suggest that it effectively transports the residual 
field produced by the decay of active regions (the 
Babcock-Leighton mechanism) towa rds the poles , 
contributing to the polarity reversal lBa bc.ocklll9fiH 
iCameron et al\ 120111: iMufioz-Jaramillo et al. 1 120 1 3l) . 
However, in the solar interior, due to the weak na¬ 
ture of this flow, its measurement is experimentally 
challenging and a complete profile mapping over the 
whole CZ has not yet been achieved by helioseismol¬ 
ogy. Since until recently the observational constraints 
on the MC profile were restricted to the surface lay¬ 
ers of the Sun, FTD models employed extrapolations 
to the solar interior based on mass conservation and 
simplified MC profiles with one or two cells per hemi¬ 
sphere. Of particular importance is the so-called re¬ 
turn flow , an equatorward directed component of the 
MC at the base of the CZ that, in the context of 
FTD, is responsible for the equatorward migration 
of the dynamo wave over the course of the cycle 
(jWang et oi.lll99lt iChoudhuri. et all 119951 ! . Even if 
more complex MC topologies are assumed, this equa¬ 
torward flo w at t h e base of the CZ retains a prepon¬ 
derant role (lHazra et thesi. 

Recently, helioseismic measurements bv lZhao et ali 
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(j2013f) mapped the MC down to a depth of 0.75 Rq 
for latitudes below 60°. Based on their measurement 
they then extrapolate the behavior of the MC all the 
way down to the base of the CZ and suggest that 
this large-scale flow is globally structured as an up¬ 
per cell that rotates counter-clockwise (in the north¬ 
ern hemisphere) and a bottom cell that rotates clock¬ 
wise. This suggested configuration implies that the 
flow direction at the base o f CZ is p olewa rd directed. 
See figures 1 and 4a of [Zhao et al\ (|2013l b However, 
when this proposed MC cell topology is used in FTD 
models, solar-like dyna mo so lutions a re n o longer pro¬ 
duced (see section 3 of lHazra et all (12014ID . Though 
these new observational developments are potentially 
transformative, we must emphasize that the actual 
measurements place little constraint on the ampli¬ 
tude and direction of the flow near the base of the 
CZ, where it matters most for regulating the cycle in 
FTD. 

Another interesting result in need of a theoretical 
interpretation is the observed variation of the surface 
MC, the amplitude of which varies in anti-phase with 
the solar cycle (Komm et al. 2009, Hathaway 2010). 
This and the development of a high-latitude counter- 
rotati ng secondary cell in the descen ding phase of the 
cycle (|Jiang et a/jl2009t lUlrichll2010l) . hints at possi¬ 
ble changes in the overall topology of the MC over 
the course of the solar cycle. Understanding how this 
process occurs requires modelling techniques that are 
dynamically consistent. 

Unlike kinematic dynamo models in which the 
flows are prescribed and their consequences stud¬ 
ied, global convective dynamo simulations self- 
consistently establish the meridional and zonal flows 
as a consequence of the interna l fluid dynamics, (e.g. 
Clune e(~atl 119991 iBruii et aU 12004 IGhizaru et all 
20101 : iKapvla et ai\ l2012f)~ In this work we use one 
such 3D MHD simulation of the solar CZ to explore 
the topological structure of the MC and its implica¬ 
tions within the context of magnetic cycles, helioseis- 
mic measurements and mean-field kinematic dynamo 
simulations. 

The data used in this study is obtained from an 
Implicit Large-Eddy Simulation (ILES) of global so- 
lar con vection prod uced by the EUL AG-MHD code 
(ISmolarkiewicz. Charbonneaul 120131 ). more specifi¬ 
cally th e millennium s i mulat ion recently presented in 
[Passos. Charbonneaul (|2014f l (henceforth PC1 4) an d 
base d on a setup describ ed in IGhizaru et id 0oioi) 
and iRacine et ali (120111 1. This simulation is per¬ 
formed on a 128x64x47 mesh in (cp,9, r) which al¬ 
lows long integration times. The numerical dissipa¬ 
tion introduced by the ILES scheme amounts to an 
implicit adaptive subgrid model which maintains nu¬ 
merical stability. In contrast to DNS, information 
on small-scale dissipation is not readily retrievable, 
but the resolved scales are dynamically consistent. 
These simulations reach a turbulent regime charac¬ 
terized by considerable dynamic range when zonal 
averages are used to define large scales and fluctu¬ 
ations around these averages are used to define small 


scales, (jRacine et oil 1201 II ). 

2. MERIDIONAL CIRCULATION STRUCTURE 

Before we investigate the dynamics of the MC in 
our simulations, we first assess the degree of resem¬ 
blance between the simula ted meridional flow and the 
solar MC. In IZhao et~cd\ (120131) the authors use he- 
lioseismic measurements carried out by SDO (HMI 
Doppler-shift data) averaged over a two year period, 
spanning from May 1, 2010 to April 30, 2012 which 
corresponds to the rising phase of sunspot cycle 24. 
For our initial assessment we will use a time window 
corresponding to the same phase of the magnetic cycle 
in our simulation. To represent solar magnetic activ¬ 
ity we built a proxy for the large-scale magnetic field 
by integrating the zonally averaged toroidal field ( B ^) 
over an extended region centered at the depth and 
latitudes where the field accumulates (see Fig. [!];>), 
following the procedure described in detail in PC14. 
This toroidal proxy is our analog of the sunspot num¬ 
ber. The large-scale magnetic cycle developing in this 
simulation shows regular polarity reversals with a pe¬ 
riod on the order of 40 years, with the large-scale 
magnetic field peaking at the base of the convection 
zone at mid/high latitudes (see figures 1 and 2 of 
PCM). 

In order to produce a latitudinal velocity diagram 
(ug) co mparable to the one presented in figure 4 of 
IZhao et aZl (1201311 . we choose a cycle from the simula¬ 
tion with similar characteristics as sunspot cycle 24, 
i.e. broadly speaking a low amplitude cycle following 
one of higher amplitude. Coincidentally, the 24 th cy¬ 
cle in our millennium simulation is a good candidate 
to perform this comparison. Figure [lp shows the pe¬ 
riod over which we averaged our signal, i.e. over the 
rising phase of the cycle_whichjs analogous to the 2 
year interval used in IZhao et ID (fMl (nearly 20% 
of the cycle duration in both cases). 

The large-scale axisymmetric latitudinal velocity 
component, ( ug ), is extracted from the simulation 
also by zonal averaging (see Fig. dJA.). Here we fol¬ 
low the convention usually used in EULAG-MHD, i.e. 
ug is defined positive towards the poles. Using ( u r ) 
and (ug), in a simila r manner to that presented in 
IGuerrero et ail (I2013T) . we also compute the stream 
function of the meridional flow which is useful for un¬ 
derstanding the cell morphology (see panels (B) and 
(C) of figure 0). 

A detailed comparison of the helioseismically- 
inferred internal meridional flow with our EULAG- 
MHD simulation is complicated by the fact that the 
large-scale magnetic cycle developing in the latter 
peaks at much higher latitude (~ 55°) than the 
sunspot butterfly diagram (~ 15°), and shows only 
hints of equatorial migration. The large-scale field 
migration toward lower latitudes in the simulation ap¬ 
pears to be inhibited by the presence of rotationally- 
aligned “banana cells” outside of a cylinder tangent 
to the equatorial base of the convection zone. The 
existence of these cells reflects the strong influence of 
rotation on turbulent convection, and is a robust fea¬ 
ture observed in all numerical simulations of rotating 
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Figure 1. (A) full disk toroidal field proxy, B. (B) Meridional representation of the toroidal field at the maximum of simulated 
cycle 24 taken at t 3 . The boxed areas mark the integration domains used to build the toroidal proxy. (C) its a zoom in on cycle 23 
an d 24 of this simul ation. The orange rectangle indicates the averaging period over which uq is computed. Compare (C) to Fig. 4 
in IZhao et al .I H2013I) . 



Figure 2. (A) large-scale ( uq ) component of the meridional flow averaged over the rising phase of simulated cycle 24 in the northern 
hemisphere. The dashed line represents the tachocline. Positive (negative) velocity values indicate poleward (equatorward) flows. 
The corresponding stream function profile that indicates the cell structure is plotted for (B) cycle minimum (at t\ on Fig.[lp) and 
(C) cycle maximum (at £ 3 ). Dashed (solid) contours in (B) and (C) represent counter-clockwise (clockwise) rotation. 


convection in the low Rossby number regime, with or the dynamical response of the simulation to the mag- 
without magnetic fields (see, e.g., Miesch & Toomre netic cycle poleward of the simulated “activity belts” 

2009, §2.2.3; also Guerrero et al. 2013). However, shows some encouraging similarities to the sun. For 
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example, Beaudoin et al. (2012) have shown that the 
pattern of cyclic rotational torsional oscillations de¬ 
veloping poleward of the mid-latitude regions of peak 
large-scale magnetic fields show the same amplitude, 
phasing, and even double-branch structure as inferred 
helioseismically (see their Fig. 5 and accompanying 
discussion). This offers hope that the meridional cir¬ 
culation dynamics may also be similar, considering 
that meridional and zonal force balances are tightly 
coupled (more on this in §3 herein). 

Figure [3] replots the mid- to high-latitude portion 
of the (ug) profile in Fig.[3]4, artificially remapped to 
the 15—60° latitude range of the iZhao et aU (1201311 
data, the latter plotted in (B) for comparison. A key 
feature of the MC extracted from the simulation is 
an equatorward flow pr esent below the largest depth 
probed bv IZhao et all (120131) . indicated by the solid 
circular arc in the left panel of figure [3] This equa¬ 
torward flow is part of a counter-clockwise meridional 
cell appearing at those latitudes and depths. Panels 
(B) and (C) of figure 2 show the different morphologi¬ 
cal configurations of the MC profile at cycle minimum 
and cycle maximum respectively. The cell where the 
return flow is located is much more pronounced at 
cyc le maximum. Moreover, as originally reported in 
iPassos et al . I (|2012H . the amplitude of this equator- 
ward flow varies in phase with the build up of the 
large-scale toroidal. This is an indication that the 
MC at this depth is highly dynamical and magnet¬ 
ically influenced. Another indication of this is the 
fact that the MC in this simulation is stronger than 
that arising in an analogous simulation using a purely 
hydrodynamic set up (no magnetic fields). 

3. MERIDIONAL CIRCULATION DYNAMICS 

The existence of the equatorward flow at mid¬ 
latitudes near the base of the convection zone is a 
robust feature of global convection simulations, and 
particularly those that include a stable radiative zone 
below the CZ (Miesch 2005). This generally occurs 
even if the MC in the bulk of the CZ is multi-celled. 
It can be attributed to two distinct but complemen¬ 
tary physical mechanisms, namely gyroscopic pump¬ 
ing and turbulent alignment of convective plumes. 
Both are essentially hydrodynamic mechanisms but 
magnetic feedbacks can alter the dynamical balances, 
inducing time variations. 

Gyroscopic pumping refers to the tendency for a 
zonal force (axial torque) to induce a meridional 
flow through th e inertia of the different ial rotation 
(|McIntvrd 120001 : iMiesch. Hindmanl 120111) . In con¬ 
vection simulations the zonal force arises principally 
from the convective Reynolds stress, which tends to 
transport angular momentum radially inward at mid¬ 
latitudes < 0) due to the Coriolis deflection 

of downflow plumes and broader upflows (Feather- 
stone & Miesch 2014). The convergence of the inward 
angular momentum flux at the base of the CZ acts to 
accelerate the zonal flow and induces an equatorward 
flow that acts to offset this acceleration. Note that 
this convergence of the convective angular momen¬ 


tum flux near the base of the CZ does not preclude 
one or more other regions of convergence closer to 
the surface. In particular, a similar convergence of in¬ 
ward angular momentum flux at the base of the near¬ 
surface shear layer due to plumes driven in the pho- 
tospheric boundary layer may account for the shallow 
return flows inferred from recent photospheric and he- 
lioseismic measurements (Hathaway 2012; Zhao et al. 
2013). In other words, the presence of an equator- 
ward return flow at ~ 0.95 Rq does not necessarily 
rule out an equatorward flow at the base of the CZ 
as well. 

Turbulent alignment refers to the tendency for 
helical downflow plumes to be diverted away from 
the vertical direction and toward the rotation axis 
(( u' r u' g ) > 0 in the northern hemisphere, opposite in 
the southern) due to the Coriolis deflection of flow 
components that are perpendicular to the rotation 
axis (Brummell et al. 1996, 1998). As these plumes 
enter the stably stratified fluid layers their vertical 
momentum is halted by the buoyancy force but they 
retain their equatorward momentum. The cumula¬ 
tive effect of many such plumes can produce a net 
equatorward circulation in the overshoot region and 
lower CZ (Miesch 2005). 

These two mechanisms for maintaining an equator- 
ward flow near the base of the CZ have been found 
previously in global convection simulations run with 
the ASH code (Miesch et al. 2000; Featherstone & 
Miesch 2014). Figure 0] demonstrates that similar 
mechanisms may be operating in our EULAG-MHD 
simulations. Shown are correlations among the fluc¬ 
tuating velocity components defined as u\ = iq — (ui), 
where Ui is the velocity component before zonal av¬ 
eraging and ( Ui) after zonal averaging. These corre¬ 
lations were computed for the northern hemisphere 
and for 2 full magnetic cycles spanning 5 polarity 
reversals (enough for illustrative purposes) over the 
integration region shown in Fig. [1)3. This region en¬ 
compasses the interface between convectively stable 
and unstable fluid layers. The (u^u'^) correlation is 
robustly negative, consistent with the maintenance 
of an equatorward circulation by the convergence of 
an inward angular momentum flux; the mechanism of 
gyroscopic pumping as discussed above. Meanwhile, 
the ( u' r u' e ) correlation is positive, indicating turbulent 
alignment. Both show variation with the phase of the 
magnetic cycle, with the inward angular momentum 
flux particularly prominent during cycle maximum. 

Though these Reynolds stresses tend to induce an 
equatorward flow near the base of the CZ, they are 
opposed at mid-latitudes by the Lorentz force. Here 
the cyclic generation of strong toroidal fields via the 
H-effect extracts energy from the rotational shear, de¬ 
celerating lower latitudes and accelerating higher lat¬ 
itudes. In the region of deceleration on the equator- 
ward edge of the toroidal bands, gyroscopic pump¬ 
ing induces a poleward meridional flow that varies in 
strength as the bands are alternately generated and 
destroyed (Fig. [2)3, C). This establishes a marked 
clockwise circulation cell at mid latitudes that ex- 
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Figure 3. (A) (ug) profile extracted from the mill ennium simulation between 40 and 75 degrees north. The dotted line represents 
the tachocline depth. (B) ug profile measured by I Zhao et all 1120131 1 between 15 and 60 degrees north. The solid black line at 
0.75 Rq represents the maximum depth probed in the measurement. The base color scheme is the same, namely warm color 
correspond to poleward flows and cold, a comparison between both images show a good resemblance in terms of the topology of this 
flow. An equatorward return flow (in blue) is clearly visible at tachocline depths in (A), a depth not yet reached by helioseismic 
measurements. 


tends into the mid-CZ. At higher and lower lati¬ 
tudes the convective angular momentum transport 
maintains a counter-clockwise cell with equatorward 
flow at the base of the CZ. This high-latitude cell 
largely disappears during cycle minima while the mid¬ 
latitude cell extends poleward (Fig. 03). 

The ( u' r u' r ) correlation in Fig. [4j which is a mea¬ 
sure of the convective overshoot, is higher when the 
toroidal field is lower and vice versa. This means 
that at this depth the penetration power of the down¬ 
flows is higher when the toroidal field is weaker. At 
cycle maximum, when the amplitude of the toroidal 
field is maximum, the overshoot layers just below the 
tachocline develops a higher “rigidity”, presumably 
due to the magnetic tension force, and the downflows 
do not penetrate as deep. As a consequence, the over¬ 
shooting downflows deposit their latitudinal momen¬ 
tum within a thinner layer immediately beneath the 
bottom of the CZ, thus driving a larger local accel¬ 
eration of the equatorward flow than at cycle mini¬ 
mum. The cyclic modulation observed in the fluctu¬ 
ating flow correlations is an important reminder that 
other forces are acting. We are currently expanding 
our analysis of these results in order to understand the 
physical mechanisms through which the large-scale 
magnetic field is modulating the turbulence. 

4. CONCLUSIONS 

Extant helioseismic measurements are now map¬ 
ping the solar meridional circulation to a depth of 
0.75 Rq. Based on 3D MHD global simulations of 
solar convection, we presented arguments that favor 
the existence of an equatorward flow below 0.75 Rq, 
near the base of the CZ at mid to high latitudes. This 
information is important for the continuous develop¬ 
ment of FTD models, and offer a challenge to future 
helioseismic measurements which will eventually re¬ 
fute or confirm this result. 


The horizontal (latitudinal) component of the MC 
extracted from this simulation, ug, is com pa red to 
the latest MC measurements of IZhao etali (|2013l) 
focusing on the latitude range poleward of the ac¬ 
tivity belts in both cases. In the simulation, solar¬ 
like features (e.g. regular large-scale magnetic cy¬ 
cle and torsional oscillations) appear at higher lat¬ 
itudes than in the Sun, and we argue that this is 
possibly due to the influence of the dynamics asso¬ 
ciated with rotationally-aligned, persistent “banana 
cells” that are typically present in this type of simu¬ 
lations. In an appropriately remapped latitude range 
and over the depths simulated (0.602 < r/R < 0.96) 
ou r ug profile sh ows the same broad characteristics 
as IZhao et aU (|2013l ) mapping presented in Fig. [3] 
We also observe evidence that the deep equatorward 
flow strength and topology is tightly connected with 
the build up of the large-scale toroidal magnetic field. 
As a possible explanation, we suggest that when this 
field component accumulates beneath the tachocline 
the associated magnetic tension force decreases the 
penetration depth of convective downflows. In turn, 
this process contributes to the enhancement of a mid 
latitude equatorward flow at the base of the CZ by 
depositing angular momentum at higher depths. 

This is a very dynamic process that changes greatly 
over the course of the solar cycle and is indicative of 
an important two-way dynamical coupling between 
field and flow. This result also suggests that spe¬ 
cial attention should be taken when interpreting the 
results provided by helioseismology, which require 
long averaging periods to extract the weak merid¬ 
ional flow from the background turbulence noise. We 
are led to believe that long periods of data integra¬ 
tion will dilute any dynamical effects impacting the 
deeper meridional flow cells. It should also be empha¬ 
sized that this tight coupling between field and flow 
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Figure 4. Time series of of the axisymmetric toroidal field, (B^) (in black), and of the Reynolds stress components (u r r 7/1), 
(u' r u'g) and (u' r u' r ), in red, green and blue, respectively. All time series are computed by integration over the boxed region of Fig. 
ITT ) in the northern hemisphere only, and are all normalized to the highest value found in the two full magnetic cycles sampled. 


should have an impact in kinematic simulations of 
flux transport dynamo models. The ubiquitous kine¬ 
matic approach is clearly missing some important as¬ 
pects of the dynamical interactions between field and 
flow. We are currently expanding the present work to 
better understand the nature of the physical mecha- 
nism(s) dynamically modulating the MC global topol¬ 
ogy in the presence of a large-scale magnetic cycle. 
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